library(RSocrata)
library(tidyverse)
library(lubridate)
library(tidygeocoder)
library(rgeoboundaries)
library(sf)
library(tmap)
library(spatstat)
source("Funciones_propias_ppp.R", encoding = "UTF-8")
La accidentalidad es uno de los problemas más comunes en el planeta provocando perdidas que pueden ir desde materiales hasta la vida misma, por lo que es de particular interés tratar de comprender como se comporta este fenómeno para evitar ser victimas de él.
La problemática en cuestión tiene la complicación de ser algo de naturaleza aleatoria, sin embargo esto no quiere decir que no exista alguna forma de controlarla pues la experiencia ha mostrado que las diferentes ciudades tienen lugares con mayor concentración de accidentes que otras ubicaciones dentro de la misma ciudad.
Debido a que se está trabajando la ocurrencia de accidentes en una ciudad de interés, la distribución Poisson y los procesos puntuales Poisson resultan ser una herramienta adecuada para tratar de explicar el fenómeno.
El estudio de accidentalidad se hará en los años 2019, 2020 y 2021 con el propósito de ver como se comporta la accidentalidad prepandemia, pandemia y postpandemia.
Una vez presentado el fenómeno que se trabajará, se selecciona la ciudad de Barranquilla para realizar el estudio.
Los datos se extrajeron de la página web GOV.CO la cual tiene diferentes bases da datos de Colombia abiertas al público. ( Accidentalidad en Barranquilla). La base de datos contiene una gran cantidad de variables, sin embargo solo se escogen aquellas de interés para el patrón puntual. La estructura de los datos es presentada a continuación:
accidentes <- read.socrata("https://www.datos.gov.co/resource/yb9r-2dsi.json") %>%
mutate(fecha = ymd(fecha_accidente)) %>%
select(12, 6:8) %>%
filter(year(fecha) > 2018 & year(fecha) < 2022)
knitr::kable(head(accidentes, 10),
col.names = c("Fecha", "Gravedad del accidentes",
"Clase de accidentes", "Sitio accidente"))
| Fecha | Gravedad del accidentes | Clase de accidentes | Sitio accidente |
|---|---|---|---|
| 2019-01-01 | Solo daños | Choque | CALLE 82 CRA 71 |
| 2019-01-01 | Con heridos | Choque | CL 30 CR 27 |
| 2019-01-01 | Con heridos | Choque | CL 17 CR 30 |
| 2019-01-01 | Solo daños | Choque | CRA 31 CALLE 68C |
| 2019-01-01 | Con heridos | Atropello | CL 99D CR 8E |
| 2019-01-01 | Solo daños | Choque | VIA 40 CR 67B |
| 2019-01-02 | Solo daños | Choque | AV CIRCUNVALAR FRENTE AL NO 70A 12 |
| 2019-01-02 | Con muertos | Choque | AV CIRCUNVALAR CR 13 |
| 2019-01-02 | Solo daños | Choque | AV CIRCUNVALAR 49 39 |
| 2019-01-02 | Solo daños | Choque | CR 53 68 216 |
Se hace necesario realizarle unos ajustes a la base de datos
considerada puesto que está no incluye las ubicaciones exactas de los
accidentes ya sea en longitud - latitud o proyección UTM, para dicho
propósito se usa la función geo() del paquete
tidygeocoder para convertir las direcciones de los
accidentes en coordenadas de longitud - latitud.
direcciones <- paste(accidentes$sitio_exacto_accidente, ", Barranquilla, Colombia", sep = "")
localizaciones <- geo(address = direcciones, method = "arcgis")
write.csv(x = localizaciones, file = "accidentes.csv", row.names = F)
Luego de realizar obtener las coordenadas en longitud - latitud, se
usan múltiples funciones del paquete sf para obtener las
respectivas ubicaciones en formato UTM obteniendose lo siguiente:
bqlla <- geoboundaries(country = "COLOMBIA", adm_lvl = 2) %>%
filter(shapeName == "DISTRITO ESPECIAL, INDUSTRIAL Y PORTUARIO DE BARR*")%>%
st_transform(crs = 3857)
coordenadas <- read.csv("accidentes.csv")
coordenadas <- cbind(accidentes, coordenadas[, 2:3]) %>%
st_as_sf(coords = c('long', 'lat')) %>%
st_set_crs(value = 4326) %>%
st_transform(crs = 3857) %>%
st_intersection(bqlla)
accidentes <- list()
for(i in 2019:2021){
accidentes[[as.character(i)]] <- filter(coordenadas, year(fecha) == i & month(fecha) > 3)
}
knitr::kable(head(accidentes[["2019"]][, c(1:4, 10)], 10))
| fecha | gravedad_accidente | clase_accidente | sitio_exacto_accidente | geometry |
|---|---|---|---|---|
| 2019-04-01 | Solo daños | Choque | CL 56 14 30 | POINT (-8327234 1227521) |
| 2019-04-01 | Con heridos | Choque | CL 30 CR 25 | POINT (-8324772 1228159) |
| 2019-04-01 | Con heridos | Choque | CR 27 CL 85 | POINT (-8330025 1229473) |
| 2019-04-01 | Solo daños | Choque | CR 38 CL 81 | POINT (-8329330 1230390) |
| 2019-04-01 | Solo daños | Choque | CALLE 31 CRA 16 | POINT (-8324126 1229674) |
| 2019-04-01 | Solo daños | Choque | CL 90 90B 42 | POINT (-8329448 1226070) |
| 2019-04-01 | Solo daños | Choque | CL 94 CR 71 | POINT (-8328610 1234554) |
| 2019-04-01 | Con heridos | Choque | AV CIRCUNVALAR CL 74 | POINT (-8331175 1231848) |
| 2019-04-01 | Con heridos | Choque | AV CIRCUNVALAR CR 13 | POINT (-8331175 1231848) |
| 2019-04-01 | Solo daños | Choque | CL 40 CR 24 | POINT (-8325400 1228241) |
tmap_mode('view')
tm_shape(bqlla)+
tm_polygons(alpha = 0.3, border.alpha = 0.7)+
tm_shape(shp = accidentes[['2019']])+
tm_dots(size = 0.01)
tmap_mode('view')
tm_shape(bqlla)+
tm_polygons(alpha = 0.3, border.alpha = 0.7)+
tm_shape(shp = accidentes[['2020']])+
tm_dots(size = 0.01)
tmap_mode('view')
tm_shape(bqlla)+
tm_polygons(alpha = 0.3, border.alpha = 0.7)+
tm_shape(shp = accidentes[['2021']])+
tm_dots(size = 0.01)
Uno de los parámetros críticos al momento de modelar un patrón puntual es la intensidad la cual puede ser constante (homogénea) o variable (inhomogénea), por lo tanto es importante iniciar el análisis con una prueba de homogeneidad de la intensidad.
# Guardando datos por anio
datos_2019 <- ppp(x = st_coordinates(accidentes[[1]])[, 1],
y = st_coordinates(accidentes[[1]])[, 2],
window = as.owin(W = bqlla))
datos_2020 <- ppp(x = st_coordinates(accidentes[[2]])[, 1],
y = st_coordinates(accidentes[[2]])[, 2],
window = as.owin(W = bqlla))
datos_2021 <- ppp(x = st_coordinates(accidentes[[3]])[, 1],
y = st_coordinates(accidentes[[3]])[, 2],
window = as.owin(W = bqlla))
Para iniciar, se divide el mapa de Barranquilla en 25 cuadrantes, si el patrón fuera homogéneo se esperaría que el número de accidentes en cada uno de las subsecciones de la ciudad contenga aproximadamente la misma cantidad de accidentes.
qc_2019 <- quadratcount(datos_2019, nx = 5, ny = 5)
plot(qc_2019, main = "Accidentes en 2019")
qc_2020 <- quadratcount(datos_2020, nx = 5, ny = 5)
plot(qc_2020, main = "Accidentes en 2020")
qc_2021 <- quadratcount(datos_2021, nx = 5, ny = 5)
plot(qc_2021, main = "Accidentes en 2021")
En todos los años se aprecia que hay una gran discrepancia entre el número de accidentes por sectores. Claramente al este de la ciudad el número de accidentes es mayor respecto al oeste lo cual es una fuerte señal de que el patrón puntual en cuestión es de naturaleza inhomogénea.
Adicional a los conteos del número de accidentes en cada uno de los sectores definidos previamente, se usa la prubea \(\chi^2\) para contrastar:
\[ \begin{cases} H_0: \text{La intensidad es constante} \\ H_1: \text{La intensidad no es constante} \end{cases} \] A continuación se presentan los resultados obtenidos
quadrat.test(qc_2019)
##
## Chi-squared test of CSR using quadrat counts
##
## data:
## X2 = 2857.1, df = 18, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 19 tiles (irregular windows)
quadrat.test(qc_2020)
##
## Chi-squared test of CSR using quadrat counts
##
## data:
## X2 = 1324.7, df = 18, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 19 tiles (irregular windows)
quadrat.test(qc_2021)
##
## Chi-squared test of CSR using quadrat counts
##
## data:
## X2 = 2094.2, df = 18, p-value < 2.2e-16
## alternative hypothesis: two.sided
##
## Quadrats: 19 tiles (irregular windows)
Puesto que en todos los casos se tienen valores - p demasiado pequeños (orden de \(10^{-16}\)) se rechaza la hipótesis de homogeneidad de intensidad en todos los casos, por lo tanto se hace necesario estimarla con algún método.
graph_ppp(datos_2019, 350, "Año 2019")
graph_ppp(datos_2020, 350, "Año 2020")
graph_ppp(datos_2021, 350, "Año 2021")